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1.  Introduction 

The  shock  response  of  solids  is  important  in  applica¬ 
tions  related  to  structural  crashworthiness,  defense  (e.g., 
projectile-target  interactions),  and  geophysics  (e.g.,  ex¬ 
plosive  mining  operations  and  hypervelocity  collisions 
of  planetary  rock  masses)  [1],  Accurate,  efficient,  stable, 
and  thermodynamically  consistent  models  for  nonlinear 
anisotropic  elasticity  are  required  for  proper  mesoscale 
modeling  of  crystalline  solids  subjected  to  impact  or  ballis¬ 
tic  loading.  Nonlinear  hyperelasticity  addresses  the  ther¬ 
modynamically  reversible  response  component  of  solids 
subjected  to  large  deformation;  classes  of  crystalline  ma¬ 
terials  of  interest  include  metals  [2],  ceramics  and  miner¬ 
als  [3],  energetic  materials  [4],  and  electronic  materials  [5], 
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This  Letter  considers  three  particular  nonlinear  thermoe¬ 
lasticity  models-each  based  on  a  different  finite  strain  ten¬ 
sor  referred  to  locally  unstressed  material  coordinates-and 
their  performance  regarding  depiction  of  shock  compres¬ 
sion  of  crystalline  solids.  The  three  strain  measures  are  the 
elastic  Green-Lagrange  strain  E  (often  simply  referred  to 
here  as  Green  strain),  Eulerian  material  strain  D,  and  loga¬ 
rithmic  material  strain  e,  all  defined  mathematically  later 
in  the  text. 

Conventional  Lagrangian  formulations  of  nonlinear 
elasticity  for  crystals  [6,7]  incorporate  the  elastic  Green 
strain  tensor.  However,  benefits  of  Eulerian  strain  ten¬ 
sors  for  isotropic  materials  were  predicted  [8]  and  demon¬ 
strated  for  cubic  crystals  under  hydrostatic  stress  [9], 
Thermal  effects  were  considered  in  [10]  for  cubic  crystals, 
and  a  theory  for  non-cubic  crystals  was  initiated  in  [11], 
A  complete  D-based  continuum  thermoelastic  theory  for 
large  deformation  of  crystals  of  arbitrary  symmetry  was 
developed  in  [12].  Analytical  solutions  for  homogeneous 
deformations  of  ideal  cubic  crystals  were  studied  over  a 
prescribed  range  of  elastic  coefficients;  stress  states  and 
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intrinsic  stability  measures  were  compared.  For  realistic 
coefficients,  Eulerian  theory  predicted  more  physically  re¬ 
alistic  and  stable  behavior  than  Lagrangian  theory  under 
large  compression  and  shear.  Analytical  solutions  for  shock 
compression  of  anisotropic  single  crystals  were  derived  for 
internal  energy  functions  quartic  in  Lagrangian  or  Eulerian 
strain  and  linear  in  entropy;  results  were  analyzed  for  non- 
metals  quartz,  sapphire,  and  diamond  in  [12]  and  metals 
aluminum,  copper,  and  magnesium  in  [  1 3].  Eulerian  theory 
was  recently  used  to  numerically  model  the  viscoplastic  re¬ 
sponse  of  aluminum  single  crystals  and  textured  polycrys¬ 
tals  in  wave  propagation  simulated  using  the  finite  differ¬ 
ence  method  [14],  wherein  Lagrangian  theory  was  found 
insufficient  for  modeling  strong/overdriven  elastic-plastic 
shocks. 

A  complete  e-based  continuum  thermoelastic  theory 
was  analogously  developed  in  [15]  and  applied  to  study 
the  shock  response  of  the  same  three  nonmetals.  The  the¬ 
ory  was  extended  to  describe  elastic-plastic  response  us¬ 
ing  a  multiplicative  decomposition  of  the  deformation  gra¬ 
dient  (formally  given  later  in  Eq.  (2)),  and  solutions  for 
plastic  shocks  (involving  slip,  twinning,  and/or  shear  frac¬ 
tures)  following  an  elastic  precursor  in  rate  independent 
solids  were  derived  [15].  Logarithmic  theory  delivered  su¬ 
perior  accuracy  to  Lagrangian  and  Eulerian  theories  for 
modeling  shocks  in  single  crystals  of  sapphire  (X-  and  Z- 
cut),  quartz  (Z-cut),  and  diamond  (X-cut)  [15],  Logarith¬ 
mic  theory  incorporating  third-order  elastic  constants  was 
also  applied  to  analytically  model  the  elastic-inelastic  re¬ 
sponse  of  isotropic  polycrystalline  titanium  diboride  ce¬ 
ramic  [16],  including  double  yield  and  effects  of  static  lat¬ 
eral  pre-stress. 

The  remainder  of  this  Letter  is  organized  as  follows. 
Kinematics  and  strain  tensors  are  formally  defined  in  Sec¬ 
tion  2,  along  with  three-dimensional  forms  of  internal 
energy  functions.  Specialization  of  the  general  theory  to 
shock  loading  and  large  volumetric  deformation  is  re¬ 
viewed  in  Section  3.  Analytical  solutions  corresponding 
to  energy  potentials  associated  with  different  strain  ten¬ 
sors  are  compared  with  each  other,  experimental  data, 
and  atomic  simulation  data  in  Section  4,  leading  to  sug¬ 
gested/preferred  potentials  for  use  in  various  situations. 
Conclusions  follow  in  Section  5.  Notation  of  continuum 
mechanics  is  used:  vectors  and  tensors  in  bold  font,  scalars 
in  italics,  and  summation  over  repeated  indices  (sub¬ 
scripts)  referred  to  Cartesian  coordinates.  As  befitting  the 
brief  style  of  a  Letter,  derivations  that  can  be  found  in  the 
cited  references  are  usually  omitted. 

2.  General  constitutive  theory 

At  a  material  element  with  reference  coordinates  X 
and  spatial  coordinates  x,  the  deformation  gradient  F  and 
volume  ratio  J  are 

F(X)  =  dx/dX, 

Fij(XK )  =  dXi/dXj  =  Sij  +  dUi/dXj ;  (1) 

J(X)  =  V/V0  =  detF; 


where  u  is  the  particle  displacement.  For  an  elastic-plastic 
material,  where  “plastic”  refers  here  to  any  thermodynam¬ 
ically  irreversible  mechanism  such  as  dislocation  glide,  de¬ 
formation  twinning,  fracture,  void  growth,  or  pore  col¬ 
lapse,  the  total  deformation  gradient  is  typically  split  into  a 
product  of  a  thermoelastic  term  (superscript  E)  and  a  plas¬ 
tic  term  (superscript  P,  and  which  can  be  further  decom¬ 
posed  into  a  product  of  deformation  mappings  associated 
with  different  physical  mechanisms)  [2-4,6,17]: 

F  =  FEFP ,  Fi]  =  FtfF^; 

J  =  JEJP  =  det  Fe  det  Fp . 

The  elastic  Green  material  strain  tensor  (i.e.,  Green- 
Lagrange  strain)  is  defined  as  [6,7] 

£=  ±(F£rF£-l),  E,J  =  ±(FEFE  -  8,j).  (3) 

Also  considered  here  are  theories  incorporating  the 
elastic  Eulerian  material  strain  tensor  [12-14] 

D  =  ^(1  -  F£_1F£_T),  Djj  =  ^(8,j  -  FlE~,FEk~^)  (4) 
and  the  elastic  material  logarithmic  strain  tensor  [15,16] 

e=  lnU£  =  -  In  (FETFE)  =  -lnC£, 

22  (5) 

e,j  =  -  \n(FETFE)ij. 

Complete  presentations  of  thermodynamic  theories  can 
be  found  in  [12-16]  and  are  too  lengthy  to  reproduce  in 
entirety  in  this  Letter.  Several  important  relations  are  listed 
next  for  reference.  Local  balance  laws  are,  in  the  absence  of 
discontinuities,  body  forces,  and  heat  conduction, 

Po  =  Pj,  V  •  or  =  px,  0  =Ja:(FF~t),  (6) 

with  p0  and  p  reference  and  spatial  mass  densities,  a  sym¬ 
metric  Cauchy  stress,  and  U  internal  energy  per  reference 
volume.  Letting  $  denote  an  internal  state  variable  linked 
to  evolution  of  microstructure  (e.g.,  defects)  and  rj  the  en¬ 
tropy  density,  forms  of  internal  energy  functions  are 

U  =  U(E,r) ,$),  U  =  U(D,  r),  §), 

U  =  U(e,  t],  |). 

Corresponding  thermoelastic  relationships  for  stress  <r 
and  temperature  6  follow  as  [12,15] 

a  =  JE~lFE(dU/dE)FET, 

or  =  JE-’FE~T(dU/dD)FE~\  (8) 

a  =JE~,FE[(dU/de)  :  (3  In  C£/9C£)]F£T; 

6  =  dO/dp,  6  =  dU/di),  6  =  dU /dp.  (9) 

At  fixed  entropy  and  internal  state  variables,  assuming 
a  stress-free  reference  configuration,  and  written  explicitly 
with  elastic  constants  up  to  fourth  order,  internal  energy 
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functions  per  unit  reference  volume  are  expressed  as  the 


Taylor  polynomials 

1  1  - 

U(E)  =  —CapEaEp  + —Ca/jyEcEpEy 

1- 

+ —CapysEaEpEyEs  +  ■  ■  ■  ,  (10) 

1  1  „ 

(J(D)  =  —  CapDaDp  +  —  CapyDaDf)Dy 

'  1  , 

+  -^CapysDaDpDyDs  +  ■  ■  ■ ,  (11) 

1  1  _ 

U  (C)  —  "T  —^E-upyGciOpQy 

'  1  _ 

"T  ~  Capys6a6p6y6§  +  ■  ■  ■  .  (12) 

Greek  indices  denote  Voigt  notation  (e.g.,  a  =  1, 


2, . . . ,  6).  Elastic  constant  tensors  of  all  orders  have  the 
same  symmetries,  though  magnitudes  of  higher-order  con¬ 
stants  differ  among  the  three  representations.  Consistent 
relationships  among  third-,  and  higher-order  elastic  con¬ 
stants  have  been  derived  elsewhere  for  various  symmetry 
classes  [9-16,18],  allowing  values  of  Lagrangian  constants 
of  potential  energy  in  (10)  to  be  converted  to  Eulerian  and 
logarithmic  analogs  in  (11)  and  (12)  without  the  need  for 
further  material  characterization  experiments.  Isentropic 
second-order  elastic  constants  Cap  are  equal  among  all 
three  theoretical  representations  in  (10)-(12). 

It  is  emphasized  that  all  theories  considered  in  the 
present  paper  are  fully  nonlinear;  no  considerations  of  nor 
comparisons  with  the  usual  linear  theory  of  elasticity  are 
included.  Sources  of  nonlinearity  are  generally  both  geo¬ 
metric,  arising  from  definitions  of  finite  elastic  strain  mea¬ 
sures  in  (3)-(5)  for  a  given  value  of  F£,  and  constitutive, 
arising  from  incorporation  of  nonzero  higher-order  elastic 
constants  (cubic  and  higher-order  terms)  in  ( 10)— ( 12).  The 
nonlinear  representations  enable  consideration  of  finite 
deformations,  of  particular  interest  large  volume  changes 
under  compressive  loading.  However,  when  elastic  dis¬ 
placement  gradient  components  and  relevant  higher-order 
constants  are  sufficiently  small  in  magnitude,  results  of  all 
three  representations  collapse  to  those  of  anisotropic  lin¬ 
ear  elasticity. 

Axial  components  of  the  three  strain  tensors  are  com¬ 
pared  for  homogeneous  spherical  and  uniaxial  elastic 
deformations  in  Fig.  1.  The  magnitude  of  the  axial  compo¬ 
nent  of  D  increases  much  more  rapidly  than  that  of  E  under 
compression,  with  e  demonstrating  trends  intermediate 
to  the  other  strains.  Internal  energy,  stress/pressure,  and 
bulk  stiffness  of  crystals  tend  to  increase  rapidly  with  com¬ 
pression  [19],  Therefore,  Eulerian  and  logarithmic  theories 
would  be  expected  to  converge  faster,  with  fewer  higher- 
order  elastic  constants  needed,  than  Lagrangian  theory,  as 
will  be  verified  later  in  this  Letter. 

3.  Shock  compression 

New  application  of  logarithmic  strain-based  theory  (e- 
based  theory)  to  shock  compression  of  metals  is  presented 
later  in  Section  4.  These  model  predictions  are  thermoe¬ 
lastic  and  strictly  applicable  only  for  very  small  volumes, 


Fig.  1.  Lagrangian  (En),  logarithmic  (en),  and  Eulerian  (Du)  strain  in 
uniform  deformation  ]  =  V /V0.  Spherical:  F  =  J1/3 1 ;  Uniaxial:  F  = 
1  +  (I  -  l)e^  ®  ei. 

such  as  in  atomic  simulations  [20,21]  or  in  the  immediate 
vicinity  of  pinned  defect  cores  [22],  wherein  plastic  defor¬ 
mation  does  not  occur.  Let  X  =  X]  denote  the  direction 
of  propagation  of  a  purely  longitudinal,  steady  thermoe¬ 
lastic  shock  wave  of  natural  velocity  £>  passing  through 
a  material  initially  unstressed  and  at  rest.  Since  the  shock 
is  modeled  as  a  jump  discontinuity  in  axial  particle  veloc¬ 
ity  v  =  3 u(X,  t)/3f,  the  strain  rate  (velocity  gradient)  in 
the  shock  front  tends  towards  infinity.  Axial  shock  stress 
P  is  defined  as  positive  in  compression  and  is  equal  to  the 
negative  of  Cauchy  stress  component  crn.  Solutions  to  this 
planar  thermoelastic  shock  problem  in  anisotropic  crys¬ 
tals  were  derived  fully  for  Lagrangian  and  Eulerian  theories 
in  [  1 2  ]  and  for  logarithmic  theory  in  [  1 5  ] ;  the  procedure  in¬ 
volves  simultaneous  solution  of  Rankine-Hugoniot  equa¬ 
tions  for  conservation  of  mass,  momentum,  and  energy: 

D  =  v/(\  —  J),  P  =  p0£>v,  U  =  “PoD2  (13) 

along  with  consideration  of  internal  energy  functions  (10), 
(11),  or  (12)  of  order  four  in  corresponding  strain  and  ex¬ 
tended  to  include  entropy  change  A  77  via  Gruneisen  tensor 

ra: 

1  1  - 

U(E,ri)  =  -CapEaEp  +  -CapyEaEpEy 

+  —CapysEaEpEyEs  +  6»0 A ^ ( 1  -  raEa),  (14) 
1  1  - 

U(D,ri)  =  -CapDaDp  + -CapyDaDpDy 

1  , 

I  ^'aPyS^cfE^pDyDg 

+  60Ar](\  -  raDa),  (15) 

1  1  , 

Lite,/))  —  —  -T  —  Capy6a6p6y 

l  6 

1  _ 

+  C  apysOa6  p6y6§  +  6*0  Af)(l  .  (16) 

The  explicit  analytical  solutions  for  shock  stress  to  order 
three  in  each  of  the  uniaxial  scalar  strain  measures  (£  = 
En,  D  =  D ii,  e  =  en)  derived  in  full  elsewhere  [12,13,15] 
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are 

P  =  -JdU/dE  =  -J[CnE+hu,E2 

+  — (Cmi -A90A)£3]  (Lagrangian),  (17) 

D 

p  =  -J~3dU/dD=  -r3[C„D+  ^CmD2 

+  -(Cmi -A<90A)D3]  (Eulerian),  (18) 

D 

p  =  -J-'dU/de=  -r’tCne+^Cme2 

+  -(Cmi -M,-Ti)e3]  (Logarithmic),  (19) 

6 

where  for  uniaxial  strain,  with  90  the  ambient  temperature 
upstream  from  the  shock: 

J  =  (1  +  2 E)1/2  =  (1  -  2D)“1/2  =  exp(e), 

„  (20) 

290A  =  Cni  +  3Cn  =  Cm  —  9Cn  =  Cm-  3Cn- 

From  (9)  and  ( 14)— ( 16),  temperature  in  the  shocked 
state  is,  for  each  theoretical  representation, 

0  =  90(  1  -  r^E)  (Lagrangian); 

6  =  80(\  -  T,D)  (Eulerian);  (21) 

6  =  (90  (1  —  Ae)  (Logarithmic). 

In  these  internal  energy-based  theories,  effects  of  en¬ 
tropy  production  on  stress  arise  from  the  rightmost  terms 
proportional  to  A  in  each  of  (17)-(19),  providing  the  in¬ 
crease  in  shock  stress  A P  on  the  Hugoniot  relative  to  that 
along  a  corresponding  isentrope.  Temperature  rise  does 
not  explicitly  enter  (17)-(  19),  but  its  effects  are  implicit  in 
the  thermodynamic  framework.  For  example,  consider  Eu¬ 
lerian  theory  with  thermoelastic  temperature  rise  A 8  = 
8  —  90  given  by  the  second  of  (21);  dissipative  stress 
increase  relative  to  temperature  rise  then  obeys  A P  = 
AD2  Ad/ 6.  Also,  noting  that  isothermal  second-order  con¬ 
stants  are  related  to  isentropic  constants  via  C0afj  =  Cap  — 
80cQrarp  with  c0  specific  heat  at  constant  strain  [6],  isen¬ 
tropic  constants  tend  to  exceed  their  isothermal  counter¬ 
parts,  leading  to  larger  P  for  adiabatic  relative  to  isothermal 
compression. 

In  laboratory  scale  specimens,  yielding  would  com¬ 
mence  in  many  pure  ductile  metals  at  small  compres¬ 
sions  at  which  effects  of  higher-order  constants  and 
differences  among  predictions  from  ( 17)— ( 19)  would  be 
negligible.  However,  nonlinear  elastic  effects  on  deviatoric 
stress  may  still  be  important  at  larger  compressions  after 
yielding,  particularly  for  lower  symmetry  materials  with 
restricted  slip  planes/directions  [13],  and  the  nonlinear¬ 
ity  in  pressure-volume  response  is  important  regardless  of 
shear  strength,  as  will  become  clear  later  in  the  context  of 
results  in  Section  4. 

When  material  strength  (i.e.,  maximum  sustained  devi¬ 
atoric  stress)  is  low  relative  to  mean  stress  and  can  safely 
be  omitted,  the  response  of  metals  to  moderate  shocks  can 
often  be  adequately  approximated  (p  ~  P)  by  an  adiabatic 
pressure-volume  ( p-J )  equation  of  state  (EOS)  derived  for 
spherical  compression,  termed  herein  the  hydrodynamic 


Table  1 

Single  crystal  and  polycrystal  properties  [13]  ( 00  =  295  K;  C ap  in  GPa;  p0 
in  g/cm3). 


Property 

Al  [100] 

Cu [100] 

Mg  [a-axis] 

Mg  [c-axis] 

Cu 

107 

166 

59.4 

61.6 

Cm 

-1080 

-1279 

-664 

-728 

Cm 

203 

715 

49 

12 

Cm 

-438 

-283 

-308 

-358 

Cmi 

25  000 

11900 

8170 

7380 

Cnn 

10  500 

2000 

1220 

893 

Cnn 

15  036 

1200 

1865 

369 

A 

2.17 

1.97 

1.52 

1.52 

Bo 

76 

137 

35.4 

35.4 

Bo 

4.42 

5.48 

3.90 

3.90 

Po 

2.70 

8.96 

1.74 

1.74 

approximation.  The  EOS  corresponding  to  each  of  (14)- 
(16)  truncated  at  order  three  in  strain  is  [12,13,15] 


-du/dj=  -B0(T1/3-;1/3) 

i  -  ^b'0(/2/3  -  i) 


(Lagrangian), 


(22) 


—  dU/dJ  =  -B0(J~7/3  —J~5/3) 

1  +  ^(b;-4)(T2/3-1) 


(Eulerian),  (23) 


—9(7/9/  =  -B0[(ln  J)/J] 
1  —  r  (Bg  —  2)  lnj 


(Logarithmic), 


(24) 


where  B0  and  B[,  are  isentropic  bulk  modulus  and  its  pres¬ 
sure  derivative  in  the  unstressed  reference  state.  Eulerian 
EOS  (23)  is  equivalent  to  the  Birch-Murnaghan  EOS  [9],  Ap¬ 
pendix  A  gives  requirements  ensuring  accuracy  of  this  hy¬ 
drodynamic  approximation  for  shock  loading. 


4.  Results  and  analysis 

Thermoelastic  and  physical  properties  for  aluminum 
(Al),  copper  (Cu),  and  magnesium  (Mg)  are  reported  in  Ta¬ 
ble  1.  Predictions  for  shock  stress  versus  volume  ratio  are 
made  using  the  analytical  solutions  of  (17),  (18),  and  (19), 
considering  only  pure  mode  directions  of  steady  shock 
propagation  (i.e.,  a  strictly  longitudinal  response,  with  no 
transverse  waves).  Shock  stress  P  normalized  by  second- 
order  isentropic  elastic  constant  Cu  is  shown  in  Fig.  2(a) 
for  Al,  Fig.  2(b)  for  Cu,  Fig.  2(c)  for  Mg  shocked  along  the  fa- 
axis]  and  Fig.  2(d)  for  Mg  shocked  along  the  [c-axis].  Cubic 
crystals  Al  and  Cu  are  shocked  along  [100]  directions. 

Higher-order  elastic  constants  for  Eulerian  and  loga¬ 
rithmic  representations  are  converted  from  measured  La¬ 
grangian  constants  using  formulae  derived  or  presented 
in  [12,13,15,18]  and  listed  in  Appendix  B.  Elastic  constants 
of  up  to  order  four  are  considered  in  results  labeled  “4th 
order”.  Results  labeled  “3rd  order”  and  “2nd  order"  are 
obtained,  respectively,  by  setting  fourth-order  and  both 
third-  and  fourth-order  elastic  constants  to  zero  in  (17)- 
(19).  Results  in  Fig.  2  compare  logarithmic  and  Eulerian 
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VIV0  VIV0 


Fig.  2.  Analytical  anisotropic  thermoelastic  solutions  for  axial  stress  in  (a)  shocked  aluminum  (Al)  single  crystal  along  [100]  compared  with  atomic 
simulation  data  [20],  (b)  shocked  copper  (Cu)  single  crystal  along  [100],  (c)  shocked  magnesium  (Mg)  single  crystal  along  a-axis,  and  (d)  shocked  Mg 
single  crystal  along  c-axis. 


Table  2 

Relative  difference  A  (%)  in  shock  stress  P  predicted  by  3rd  (Lag3,  Log3,  Eul3)  and  2nd  (Lag2,  Log2,  Eul2)  order  theories:  Al  and  Cu. 


v/v0 

Aluminum  [100] 

Copper [100] 

Lag3 

Log3 

Eul3 

Lag2 

Log2 

Eul2 

Lag3 

Log3 

Eul3 

Lag2 

Log2 

Eul2 

0.96 

-5.0 

-3.8 

-3.1 

-23.1 

-11.9 

+  1.1 

-1.6 

-0.2 

-0.4 

-15.8 

-3.7 

+9.3 

0.92 

-15.7 

-13.9 

-14.1 

-48.2 

-29.8 

-5.0 

-5.4 

-0.8 

-2.1 

-31.6 

-7.8 

+  19.9 

0.88 

-28.1 

-28.2 

-34.5 

-71.2 

-51.0 

-19.8 

-10.4 

-1.8 

-6.2 

-46.5 

-12.4 

+31.5 

Table  3 

Relative  difference  A  in  shock  stress  P  predicted  by  3rd  (Lag3,  Log3,  Eul3)  and  2nd  (Lag2,  Log2,  Eul2)  order  theories: 

Mg. 

V/V0 

Magnesium 

i  [a-axis] 

Magnesium 

[c-axis] 

Lag3 

Log3 

Eul3 

Lag2 

Log2 

Eul2 

Lag3 

Log3 

Eul3 

Lag2 

Log2 

Eul2 

0.96 

-2.9 

-0.8 

-0.6 

-22.8 

-10.9 

1.2 

-2.5 

-0.2 

-0.4 

-23.4 

-11.0 

-0.0 

0.92 

-9.2 

-3.0 

-2.9 

-44.7 

-22.7 

1.0 

-7.9 

-0.6 

-2.0 

-45.1 

-22.5 

-1.2 

0.88 

-16.8 

-6.3 

-7.6 

-64.3 

-35.0 

-1.4 

-14.5 

-1.2 

-5.2 

-64.2 

-33.0 

-3.8 

theories;  different  plots  comparing  Lagrangian  and  Eule- 
rian  theories  can  be  found  in  [13],  All  longitudinal  higher- 
order  elastic  constants  (i.e.,  ail  third-  and  fourth-order  con¬ 
stants)  are  smaller  in  magnitude  for  Eulerian  and  logarith¬ 
mic  theory  than  Lagrangian  theory  for  these  metals,  as  can 
be  verified  via  examination  of  values  in  Table  f . 

Stress  predictions  of  2nd  and  3rd  order  models  are  usu¬ 
ally  closer  to  those  of  4th  order  theory  for  Eulerian  and  log¬ 
arithmic  theory  than  for  Lagrangian  theory.  Tables  2  and  3 
list  relative  differences  (%)  of  2nd  and  3rd  order  predictions 
compared  to  4th  order  predictions,  computed  as 

A  =  2  •  (2nd  or  3rd  order  result  —  4th  order  result)/ 


(2nd  or  3rd  order  result  +  4th  order  result).  (25) 

For  each  crystal  type,  such  differences  are  almost  al¬ 
ways  smaller  in  magnitude  for  Eulerian  theory  (Eul2,  Eul3) 
and  logarithmic  theory  (Log2,  Log3)  than  for  Lagrangian 
theory  (Lag2,  Lag3)  at  a  given  volume  ratio  and  order  of 
approximation,  i.e.,  order  of  Taylor  polynomial  used  in  in¬ 
ternal  energy  functions  (14)-(16).  If  4th  order  results  are 
viewed  as  exact,  then  A  can  be  interpreted  as  a  measure 
of  error  of  the  lower  order  predictions;  however,  the  term 
“error”  is  avoided  in  the  present  context  since  even  the 
4th  order  predictions  may  be  subject  to  inaccuracy  relative 
to  “exact”  polynomial  representations  containing  up  to  an 
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infinite  number  of  terms  or  relative  to  experimental  data 
used  for  model  validation. 

These  observations  imply  a  faster  converging  series  in 
(19)  or  (18)  than  (17)  as  the  number  of  higher-order  con¬ 
stants  is  increased  and  presumably  greater  accuracy  of  Eu- 
lerian  or  logarithmic  theory  than  Lagrangian  theory  when 
the  same  number  of  constants  (i.e.,  the  same  order  of  Taylor 
polynomial)  is  used  in  each  representation.  Eulerian  and 
logarithmic  theories  trend  towards  comparable  accuracy, 
with  one  theory  or  the  other  possibly  more  accurate  for  a 
given  metal  or  order  of  approximation.  For  Al,  as  shown 
in  Fig.  2(a),  the  analytical  solution  incorporating  4th  order 
Eulerian  theory  best  matches  atomic  predictions  [20],  For 
example,  at  a  compression  ratio  of  J  =  V /V0  =  0.923, 
the  difference  in  4th  order  Eulerian  theory  versus  atomic 
simulation  is  —1.4%,  compared  to  —7.4%  for  4th  order  La¬ 
grangian  theory  and  —2.2%  for  4th  order  logarithmic  the¬ 
ory. 

Comparison  with  results  of  empirical  potential-based 
molecular  dynamics  (MD)  calculations  [20]  may  be  viewed 
with  some  reservation  since  MD  results  are  subject  to  their 
own  possible  sources  of  inaccuracy,  e.g.,  uncertainty  in  as¬ 
sumed  functional  forms  and  parameters  of  the  potentials 
and  precision  of  numerical  solutions.  However,  results  of 
the  Winey-Kubota-Gupta  (WKG)  potential  considered  in 
Fig.  2(a),  first  developed  in  [2 1  ],  were  reported  in  [20]  to  be 
more  suitable  than  those  of  several  other  embedded  atom- 
type  potentials  for  modeling  the  uniaxial  shock  response  of 
defect-free  aluminum  single  crystals  and  in  [21]  more  ac¬ 
curate  for  modeling  thermoelasticity  in  Al  over  a  wide  tem¬ 
perature  range.  Furthermore,  later  in  this  paper,  predic¬ 
tions  of  the  nonlinear  continuum  models  are  compared  di¬ 
rectly  with  data  from  shock  experiments  and  density  func¬ 
tional  theory  (DFT)  calculations  (e.g.,  labeled  as  “ab  ini¬ 
tio",  “first  principles",  “electronic  structure”  or  “quantum- 
MD”),  the  latter  typically  viewed  as  less  empirical  than 
classical  MD  calculations,  though  again  various  assump¬ 
tions  and  sources  of  inaccuracy  are  involved  in  solution  of 
the  governing  quantum-mechanical  relations  [7], 

Predicted  shock  stress  P  need  not  always  increase 
(or  decrease)  as  the  order  of  approximation  used  in  any 
of  (17)-(19)  is  increased.  Under  uniaxial  compression, 
finite  strain  measures  E,  D,  and  e  are  all  negative  in 
sign,  leading  to  compressive  (positive)  contributions  to  P 
from  positive  second-order  constants,  tensile  (negative) 
contributions  to  P  from  positive  third-order  constants,  and 
compressive  (positive)  contributions  to  P  from  positive 
fourth-order  constants.  Signs  of  these  contributions  are 
reversed  when  signs  of  elastic  constants  are  negative.  For 
materials  and  properties  of  Table  1,  all  second-  and  fourth- 
order  elastic  constants  are  positive;  third-order  constants 
are  all  negative  for  Lagrangian  and  logarithmic  theory  and 
are  all  positive  for  Eulerian  theory.  As  a  result,  predictions 
for  logarithmic  theory  in  Fig.  2  increase  monotonically  with 
increasing  order,  while  predictions  of  Eulerian  theory  do 
not,  with  the  3rd  order  Eulerian  prediction  of  P  smaller 
than  that  of  the  2nd  order  Eulerian  prediction.  As  a 
general  implication  with  regard  to  arbitrary  materials,  both 
the  sign  and  magnitude  of  higher-order  elastic  constants 
should  be  considered  when  interpreting  their  importance 
in  the  context  of  stress  predictions. 


v/v0 


Fig.  3.  Analytical  solutions  for  polycrystalline  (a)  aluminum,  (b)  copper, 
and  (c)  magnesium  in  the  hydrodynamic  limit  and  experimental  data  [23]. 

Removing  prior  stated  restrictions  to  defect-free  or 
small  volumes  of  single  crystals,  predictions  of  each  third- 
order  EOS  in  (22)-(24)  are  compared  with  polycrystalline 
shock  compression  data  compiled  in  [23]  for  Al  in  Fig.  3(a), 
Cu  in  Fig.  3(b),  and  Mg  in  Fig.  3(c).  Compressibility  proper¬ 
ties  of  Table  1  used  in  these  predictions  are  obtained  from 
ultrasonic  experiments  [24]  and  are  not  fit  to  the  shock 
data. 

Individual  references  from  which  experimental  data  are 
complied  include  [25-27]  for  Al,  [25,28]  forCu.and  [26]  for 
Mg.  Experimental  methods,  as  discussed  in  detail  in  [25— 
29],  involve  measurements  of  shock  wave  velocity  and  in 
some  cases  free  surface  velocity  for  a  specimen  subjected 
to  a  planar  shock  generated  by  detonation  in  a  high  explo¬ 
sive  system.  Maximum  uncertainties  in  volume  change  at 
fixed  Hugoniot  stresses  as  reported  in  these  works  range 
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from  ±0.5%  to  ±2%  [25-29].  Error  bars  on  p/B0  at  fixed  V 
are  reasonably  translated  to  ±3%  in  Fig.  3. 

Relevance  of  the  hydrodynamic  assumptions  detailed  in 
Appendix  A  follows  from  the  relatively  low  dynamic  yield 
strength  Y  of  these  high  purity  metals,  with  values  on  the 
order  of  tens  to  several  hundred  MPa  for  dislocation  glide 
in  A1  [14,18]  and  Cu  [18,29]  and  for  slip  or  twinning  in 
Mg  [30-32].  At  very  high  rates,  although  Y  may  increase 
with  pressure  or  decrease  with  thermal  softening  and 
shear  localization,  it  follows  that  for  V /VQ  <  0.95,  p/B0  > 
0.05  in  Fig.  3  and  Y/ B0  ~  O(10-3),  leading  to  Y/p  < 
0.02  <?C  1  for  each  metal  in  these  regimes  of  interest.  The 
metallic  polycrystals  are  effectively  isotropic  and  plasti¬ 
cally  incompressible  (i.e.,  non-porous).  Hydrodynamic  ap¬ 
proximations  were  also  found  accurate  for  shocked  A1  and 
Cu  in  other  works  [18,29,33]. 

In  each  part  of  Fig.  3,  Eulerian  theory  is  most  accurate, 
followed  by  logarithmic  theory  and  then  Lagrangian 
theory,  with  the  latter  much  too  compliant.  Because 
Eulerian  theory  appears  superior  for  modeling  spherical 
compression  and  at  least  equally  valid  as  logarithmic 
theory  for  modeling  uniaxial  strain  compression,  Eulerian 
theory  is  herein  deemed  preferable  overall  for  representing 
the  shock  response  of  ductile  metals. 

As  demonstrated  in  Appendix  C,  Eulerian  theory  also 
compares  favorably  with  DFT  results  for  static  isothermal 
compression  of  Al,  Cu,  and  Mg,  especially  with  regard  to 
relative  accuracy  of  lower-order  EOS  predictions.  Future 
numerical  studies  should  thus  consider  applying  Eulerian 
theory  in  finite  element  or  finite  difference  simulations 
of  the  response  of  ductile  metals,  offering  potential 
improvement  over  prior  implementations  incorporating 
Lagrangian  theory  [34-36]. 

On  the  other  hand,  recent  work  [15]  found  that  loga¬ 
rithmic  theory,  while  more  computationally  cumbersome 
than  Eulerian  theory,  more  accurately  captures  the  shock 
response  of  ceramic  crystals  sapphire  (ct'-AbOs),  diamond 
(C),  and  quartz  (a-Si02 ),  requiring  fewer  higher-order  elas¬ 
tic  constants  to  attain  close  agreement  with  experimental 
plate  impact  data.  As  shown  in  Table  4,  these  nonmetals 
have  a  large  ratio  of  effective  shear  (G0)  to  bulk  modulus 
compared  to  metals,  characteristic  of  the  trend  that  G0/B0 
tends  to  increase  with  brittleness  in  pure  substances  [37— 
39],  This  correlation  between  plastic  malleability  and  de¬ 
creasing  G0/B0  =  3(1  —  2v)/(2  +  2v),  first  proposed 
by  Pugh  [37],  is  explained  in  terms  of  atomic  structure 
and  bonding  in  [38]  and  verified  in  the  context  of  inden¬ 
tation  hardness  in  [39].  Generally,  resistance  to  bending  of 
atomic  bonds  is  reflected  by  the  shear  modulus,  and  re¬ 
sistance  to  stretching  or  compression  of  atomic  bonds  by 
the  bulk  modulus  [38],  Second-order  elastic  constants  can 
be  related  to  harmonic  interatomic  forces,  whereas  higher- 
order  constants  depend  on  anharmonicity  [40].  Future  nu¬ 
merical  studies  of  deformation  of  rock  such  as  granite  [41] 
(silicate  minerals)  should  benefit  from  use  of  logarithmic 
rather  than  Lagrangian  nonlinear  elasticity. 

Since  the  nonmetals  considered  in  Table  4  tend  to  have 
a  rather  large  Hugoniot  Elastic  Limit  (HEL)  (e.g.,  not  oc¬ 
curring  until  up  to  »10%  compression)  and  may  demon¬ 
strate  significant  shear  strength  when  shocked  above  their 


Table  4 


Ratio  of  ambient  Voigt-averaged  shear  modulus  to  bulk  modulus  [13,15]. 


Cu 

Al 

Mg 

cy-A^C^ 

C 

Qf-Si02 

Gq/Bq 

0.34 

0.35 

0.53 

0.65 

1.22 

1.27 

HEL  [12],  use  of  a  pressure-volume  EOS  alone  is  insuffi¬ 
cient  for  calculation  of  axial  shock  stress.  However,  suit¬ 
ability  of  logarithmic  EOS  (24)  for  representing  the  hydro¬ 
static  pressure-volume  response  has  been  verified  for  alu¬ 
mina,  diamond,  and  quartz  in  [15],  Table  5  summarizes 
loading  conditions,  materials,  and  model  performance 
studied  in  the  present  work  (*)  as  well  as  results  in  [12, 
13,15],  In  summary,  Eulerian  nonlinear  theory  is  recom¬ 
mended  for  modeling  the  response  of  ductile  metals,  while 
logarithmic  theory  is  recommended  for  modeling  ceramics 
and  minerals  with  larger  ratios  of  shear  to  bulk  stiffness. 

New  contributions  of  the  present  paper  are  clarified 
as  follows.  Complete  thermomechanical  theories  incorpo¬ 
rating  Eulerian  and  logarithmic  strain  measures  were  for¬ 
malized  in  prior  work  [12—15];  the  present  paper  contains 
new  results  and  comparisons  of  applications  of  such  theo¬ 
ries  for  three  metals  (Al,  Cu,  and  Mg)  subjected  to  uniaxial 
shock  compression  (single  crystals)  and  volumetric  com¬ 
pression  (polycrystals).  In  particular,  predictions  of  Eule¬ 
rian  and  logarithmic  theories  for  these  metals  had  not  been 
compared  previously  with  respect  to  order  of  representa¬ 
tion,  experimental  data,  or  atomic  simulation  data  (neither 
MD  nor  DFT).  Such  new  comparisons,  when  viewed  in  con¬ 
junction  with  elastic  constants  listed  in  Table  4  and  prior 
results  listed  in  Table  5,  have  enabled  discovery  of  a  new 
inverse  correlation  between  Pugh’s  ductility  ratio  [27]  and 
suitability  of  Eulerian  versus  logarithmic  theory  for  mod¬ 
eling  shock  compression  processes,  with  the  former  pre¬ 
ferred  for  lower  values  of  the  ratio  Go/B0. 

5.  Conclusions 

Constitutive  theories  and  analytical  solutions  for  ther¬ 
moelastic  compression  relevant  to  shock  loading  (extreme 
strain  rates  and  pressures)  have  been  presented.  Theories 
consider  Taylor-type  polynomial  representations  of  inter¬ 
nal  energy  density  in  one  of  three  finite  strain  measures 
referred  to  unstressed  material  coordinates:  Green  strain, 
Eulerian  strain,  or  logarithmic  strain.  The  present  results, 
newly  examined  in  conjunction  with  prior  work,  imply 
that  Eulerian  theory  should  offer  the  greatest  accuracy  (at 
fixed  order  of  Taylor  polynomial)  for  modeling  the  shock 
response  of  typical  ductile  metals.  In  contrast,  logarithmic 
theory  is  implied  as  most  accurate  for  representation  of 
the  shock  response  of  more  brittle  ceramics  and  minerals 
characterized  by  a  relatively  larger  ratio  of  effective  elastic 
shear  modulus  to  bulk  modulus  in  the  initial  state. 

Appendix  A.  Hydrodynamic  approximation 

For  the  hydrodynamic  approximation  invoked  implic¬ 
itly  in  application  of  EOS  calculations  to  shock  response, 
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Table  5 

Summary  of  present  (*)  and  prior  research  results:  finite  strain  model  evaluations. 


Loading  Protocol 

Material 

Recommended  theory 

Remarks 

Ref. 

Hydrostatic  compression 

Ideal  cubic,  Bq  =  4 

Eulerian 

Eulerian  more  accurate  p-V  curves 

[12] 

Uniaxial  compression 

Ideal  cubic,  Bq  =  4 

Eulerian 

Eulerian  more  accurate  and  stable 

[12] 

Simple  shear 

Ideal  cubic,  Bq  =  4 

Eulerian 

Eulerian  more  accurate  and  stable 

[12] 

Shock  compression 

Sapphire  (AI2O3) 

Logarithmic 

Logarithmic  more  accurate  overall 

[15] 

Shock  compression 

Diamond  (C) 

Logarithmic 

Logarithmic  more  accurate  overall 

[15] 

Shock  compression 

a-Quartz  (Si02) 

Logarithmic 

Logarithmic  more  accurate  (Z-cut) 

[12,15] 

Shock  compression 

Aluminum  (Al) 

Eulerian 

Eulerian  best  fit  to  atomic  data 

[13,*] 

Shock  compression 

Copper (Cu) 

Eulerian 

Eulerian  faster  convergence 

[13,*] 

Shock  compression 

Magnesium  (Mg) 

Eulerian 

Eulerian  faster  convergence 

[13,*] 

Shock  compression 

Al  polycrystal 

Eulerian 

Eulerian  more  accurate 

[*] 

Shock  compression 

Cu  polycrystal 

Eulerian 

Eulerian  more  accurate 

[*] 

Shock  compression 

Mg  polycrystal 

Eulerian 

Eulerian  more  accurate 

[*] 

deviatoric  stresses  are  assumed  negligible  relative  to  hy¬ 
drostatic  stress: 

a  «  -pi  =>.  o-n  sa  o22  %  (r33  «  -p 

=>•  P  =  -on  ^  P-  (A.1) 

Contributions  from  stress  to  internal  energy  in  the  third 
of  (6)  and  (13)  arise  only  from  pressure: 

U  =Jcr  :  Vu  =  —JpV  ■  v  =  —pj, 
p  =  -dU/dJ  =  -VodU/dV,  (A.2) 

where  plastic  volume  change  has  been  assumed  negligible 
in  the  second  of  (A.2)  such  that  J  =  JE.  Total  deformation 
need  not  be  spherical  in  (A.2),  but  only  volume  changes 
perform  mechanical  work.  Use  of  (A.l),  (A.2)  and  (B.3)  in 
(14)— ( 16)  yields  the  EOS  representations  in  (22)-(24).  Use 
of  (9)  in  (14)-(16)  with  3T0  =  trr  yields  the  temperature 
rise  resulting  from  a  volume  change: 


•  The  material  possesses  cubic  or  isotropic  elastic  sym¬ 
metry  for  (B.3)  to  hold,  with  isotropy  applicable  to  poly¬ 
crystalline  metallic  samples  with  many  randomly  ori¬ 
ented  grains; 

•  Contributions  of  fourth-  and  higher-order  in  elastic 
strain  to  U  (and  thus  orders  three  and  higher  to  P)  can 
be  omitted,  which  in  turn  leads  to  omissions  of  higher- 
order  derivatives  such  as  B'0'  and  A  in  each  EOS  in  (22)- 
(24).  Adequacy  of  this  truncation  can  only  be  judged  via 
accuracy  of  predictions  of  P  orp,  e.g.,  with  experimental 
data. 

Only  the  first  of  the  above  five  points  is  mandatory  for 
any  general  “hydrodynamic  theory”  that  neglects  material 
rigidity  [25]. 

Appendix  B.  Higher-order  elastic  constants 


e  =  e0 

e  =  eQ 
9  =  e0[\ 


i  +  -r0(i -j2/3) 


(Lagrangian); 


1  +  -r0(j~2/3  -  i) 


(Eulerian); 
—  /~o  Injj  (Logarithmic). 


(A.3) 


The  following  modeling  assumptions  are  inherent  in  the 
above  derivations  and  the  current  application  of  (22)-(24) 
to  describe  planar  shock  compression  processes  in  ductile 
metals: 


•  Deviatoric  stress  components  are  negligible  relative  to 
pressure;  since  such  components  are  on  the  order  of 
dynamic  yield  stress  Y  for  metals  [29],  this  translates 
to  Y/p  1; 

•  Entropy  contributions  to  U  can  be  adequately  repre¬ 
sented  by  rightmost  linear  terms  in  ( 14)— ( 16);  this  re¬ 
quires  that  magnitudes  of  entropy  increases  across  the 
shock  are  modest,  which  is  in  accord  with  the  first  ap¬ 
proximation  above  since  plastic  work  contributions  to 
entropy  rise  are  correspondingly  small  when  deviatoric 
stresses  are  sufficiently  small; 

•  Other  possible  sources  of  internal  energy  change  such 
as  resulting  from  heat  conduction,  defect  generation 
(£),  and  plastic  volume  change  (/p)  are  negligible;  the 
first  is  in  accord  with  usual  adiabatic  assumptions 
inherent  in  (6)  and  (13),  the  last  in  accord  with  usual 
kinematic  descriptions  that  plastic  deformation  due  to 
slip  and  twinning  are  isochoric; 


Relationships  among  higher-order  elastic  constants  en¬ 
tering  Lagrangian,  Eulerian,  and  logarithmic  representa¬ 
tions  have  been  derived  in  several  ways,  including  match¬ 
ing  terms  of  strain  energy  potentials  of  equivalent  pow¬ 
ers  [11,18],  matching  of  tangent  elastic  moduli  [12],  and 
use  of  known  formulae  for  derivatives  of  proper  symmetric 
fourth-order  tensor  functions  [42].  Consistent  results  from 
all  methods  of  derivation  for  third-order  longitudinal  con¬ 
stants  for  a  stress-free  reference  configuration  are  [11-13, 
18,42] 

fill  =  Cm  —  12Cn  =  Cm  —  6Cm  (B.l) 

For  fourth-order  longitudinal  elastic  constants  [13,18], 
Cnn  =  fm  +  18Cm  +318Cn 

=  fun  —  12Cm  —  28Cm  (B.2) 

The  first  pressure  derivative  of  the  bulk  modulus  at  the 
reference  state,  for  crystals  of  cubic  or  isotropic  symmetry, 
with  3B0  =  C,i  +  2C12  obeys  [6,12,13,15,16] 

B0  =  “  ^gflll  +  -C112  +  -Cl23^  j  Bo 

=  —  ^  g^111  9^123)  f  ^ 

=  —  Cm  +  -  C112  +  -  Ci23^  j B0  +  2.  (B.3) 
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Fig.  4.  Analytical  EOS  solutions  and  DFT  data  [43-45]  for  pure  single 
crystalline  (a)  aluminum,  (b)  copper  (c)  magnesium. 

Appendix  C.  Comparison  with  static  DFT  pressure- 
volume  data 

Predictions  of  each  EOS  in  (22)-(24)  are  compared 
with  isothermal  compression  data  obtained  from  DFT 
calculations  reported  in  [43-45]  for  pressures  up  to 
100  GPa.  Results  for  A1  and  Cu  (FCC  phases)  correspond  to 
room  temperature  (293  1<)  [43,44];  for  Mg  (F1CP  phase)  to 
0  K  [45].  For  proper  comparison,  isothermal  values  of  bulk 
modulus  and  its  pressure  derivative  as  reported  from  the 
DFT  works  for  each  type  of  perfect  single  crystal  are  used 
in  each  EOS  rather  than  isentropic  polycrystalline  values 
given  previously  in  Table  1.  Specifically,  continuum  EOS 
predictions  labeled  “3rd  order"  invoke  values  of  B0  and  B[, 
shown  in  the  inset  of  each  part  of  Fig.  4.  Continuum  EOS 
predictions  labeled  "2nd  order”  correspond  to  omission,  in 


each  EOS,  of  terms  involving  third-order  elastic  constants 
in  (B.3),  implying  values  of  zero,  two,  and  four  for  B[,  in  (22), 
(24),  and  (23),  respectively,  noting,  however,  that  (B.3) 
does  not  strictly  apply  for  single  crystal  Mg  with  hexagonal 
symmetry. 

As  shown  in  Fig.  4,  3rd  order  Eulerian  and  logarithmic 
theories  are  of  comparable  accuracy  relative  to  DFT  data 
for  Al,  Cu,  and  Mg,  with  3rd  order  Lagrangian  theory 
less  accurate  (e.g„  too  compliant)  at  volume  ratios  J 
less  than  «0.85.  Second-order  logarithmic  theory  is  too 
compliant  for  each  material  at  volume  ratios  less  than 
«0.90,  and  2nd  order  Lagrangian  theory  is  even  more 
so,  with  the  lowest  relative  accuracy  among  the  three 
theories  for  the  present  application.  Flowever,  2nd  order 
Eulerian  theory  is  of  comparable  accuracy  to  3rd  order 
Eulerian  and  to  3rd  order  logarithmic  theory  for  Al  and 
Mg  (Fig.  4(a)  and  (c))  and  is  more  accurate  than  the  two 
other  2nd  order  continuum  elasticity  predictions  for  Cu 
(Fig.  4(b)).  In  summary,  Fig.  4  further  supports  the  present 
recommendation  for  preferred  use  of  Eulerian  elasticity  for 
modeling  the  high  pressure  response  of  ductile  metals. 
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